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Abstract 



We investigate fluids of dipolar hard particles by a certain variant of 
density-functional theory. The proper treatment of the long range of the dipo- 
lar interactions yields a contribution to the free energy which favors ferromag- 
netic order. This corrects previous theoretical analyses. We determine phase 
diagrams for dipolar ellipsoids and spherocylinders as a function of the aspect 
ratio of the particles and their dipole moment. In the nonpolar limit the re- 
sults for the phase boundary between the isotropic and nematic phase agree 
well with simulation data. Adding a longitudinal dipole moment favors the 
nematic phase. For oblate or slightly elongated particles we find a ferromag- 
netic liquid phase, which has also been detected in computer simulations of 
fluids consisting of spherical dipolar particles. The detailed structure of the 
phase diagram and its evolution upon changing the aspect ratio are discussed 
in detail. 

PACS numbers: 64.70.Md, 61.30.By, 77.80.-e 
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I. INTRODUCTION 



There are two basic molecular properties which can cause long-ranged orientational order 
in fluids. First, as has already been shown by Onsager particles of sufficiently anisotropic 
shape, e.g. long rods or flat discs, form a nematic phase at high densities. This phase transi- 
tion can be induced by purely steric interactions This has been confirmed by computer 
simulations of hard ellipsoids [||, spherocylinders and cut spheres [H],nj] ; which have 

become standard models of liquid crystals. Some of these systems exhibit further transitions 
to a smectic or a columnar phase. Second, there is growing evidence that a ferromagnet- 
ically ordered nematic phase can be stabilized by dipolar interactions between spherical 
particles, i.e. in the absence of anisotropic steric interactions. This phase has been observed 
in Monte Carlo simulations of dipolar soft Il2|,[l3j and hard [|l4|-p^] spheres as well as in 



Stockmayer fluids fll7[ . It has also been analyzed by density-functional theory [|18| -|20||. Due 
to the long range of the dipolar interactions in this phase the equilibrium configuration 
exhibits a spatially inhomogeneous magnetization similar to the domain formation in 
solid ferromagnets. 

Molecules typically possess both a shape anisotropy and a permanent dipole moment. 
Therefore it is interesting to analyze the relative importance of these two properties with 
respect to the formation of orient at ionally ordered phases and the crossover from a ferro- 
magnetic to a purely nematic phase. To this end in the present work we study the models of 
dipolar hard ellipsoids and spherocylinders which cover the models of dipolar hard spheres 
as well as of nonpolar elongated or oblate hard particles as limiting cases. These models have 
already been examined by Onsager's virial expansion |p2 |, integral equation theories |23| . 



and different kinds of density-functional theory E4,25|. However, some of these approaches 



22,23|j25|| suffer from an incorrect treatment of the long-ranged dipolar forces and therefore 



fail to predict a ferromagnetic phase. In simulations of dipolar hard ellipsoids and sphe- 
rocylinders |^],[^,|rj]] and the dipolar Gay-Berne model |29| this phase has also not yet been 
found, probably because the simulations were restricted to large elongations of the particles 
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which tend to destabilize the ferromagnetic order. Terentjev et al. [^] find theoretical in- 
dications that a ferromagnetic phase might form more readily in dipolar liquid- crystalline 
polymers. 

In the present work we examine the models mentioned above by an alternative density- 
functional theory which is a generalization of the theory applied to Stockmayer fluids in 
Refs. §}WfZl\. 



II. DENSITY-FUNCTIONAL THEORY 

As motivated in the introduction, we consider fluids consisting of hard particles which 
have a symmetry axis and carry in their center a pointlike permanent dipole moment aligned 
with this axis. The orientation of these uniaxial particles with respect to spatially fixed 
coordinates is described by two angles (9, <fi) = lj. The interaction pair potential u>(ri 2 , lj, lj') 
is the sum of the hard core potential and the dipolar potential: w = Whc + Wdi P - The former 
is given by 

oo, r 12 < a(u 12 ,u,uj') 
w hc {ri2,uj,uj) = < (1) 

I 0, otherwise 

where r 12 = r — r' is the center-to-center distance vector between the two particles at r and 
r', respectively, and a{ui 2 , to, lj') is the distance of closest approach for given orientations lj, 
lj', and lj\ 2 of the particle axes and the vector r 12 , respectively. The dipolar potential has 
the form 

2 

TfL 

w dip {r l2 , lj, lj') = -— [3(m(o;) • ri 2 )(m(^') • r 12 ) - m(oj) ■ m(a/)] (2) 

2 

m _ , , 



,3 



W\LJ 12 ,LJ,LJ 



In Eq. (0) m(u) is the dipole vector and m is its absolute value. The hats denote unit 
vectors. 
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A. Reference system of hard particles 



We first analyze the reference system of the corresponding nonpolar hard core fluid. Its 
free energy as a functional of the number density p(r, uj) of particles at r with orientation uj 
can be written as 

Fref = F id + F ex (3) 

where the ideal gas part is given by (here and in the following we use the notations J d 3 r = 
J v d 3 r and f duo = f s duj where V denotes the volume of the sample and S 2 the unit sphere) 

(3F ld = J d 3 rdLup(r,u){\n(47rp(r,uj)\ 3 ) - 1], (4) 

A being the thermal wavelength. The excess part of the free energy is related to the direct 
correlation function c(t,t',uj,uj'; [p]): 

~ r W, I , = c(r,r>,u/; [p]). (5) 

k B T op(r, u)dp{r, uj') 

Equation (g) can be integrated twice along a linear path in density space starting from a 
zero density state yielding pi 



(3F ex = - j d 3 rdud 3 r'du' J d\(l - A)c(r, r', uj, uj'] [Xp\)p(r, uj)p(r', uj'). (6) 

Due to the absence of exact results in order to proceed one now needs an approximation 
for the direct correlation function. An educated guess, which renders a computationally 
simple approach but nevertheless yields reliable results for the isotropic-nematic transition 



of nonpolar hard particles, is given by the decoupling approximation introduced by Pynn |32 
which assumes that this anisotropic function can be obtained from a dimensionless function 
co(x;rj) of a single variable by a suitable anisotropic rescaling with the distance of closest 
approach: 

c(r,r',uj,uj';[p]) w c (r 12 /cr(uj 1 2,uJ,uj');rj). (7) 

Here 77 = pv^ ' is the packing fraction of the particles with an individual volume and 
p = -^y j d 3 rdup(r, uj) is the mean number density. In this work we confine ourselves to 
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spatially homogeneous phases so that p(r,u) = pa(u) with the normalized orientational 
distribution a(u))\ J dua(u) = 1. This implies that we do not consider smectic or solid 
phases and that in the case of ferromagnetic order the sample shape is taken to be needlelike 
which suppresses the formation of domains PTl . With these assumptions and approximations 



Eq. (^|) reduces to 



V 



P fo(v) / dujduj'duji<2a{uj)a{uj')(j (uj\2, uj, uj'). (8) 



The function Co enters only in the integrated form 



poo i>l 

f (rj) = - dxx 2 / d\(l - A)c (x; Xr}). (9) 
Jo Jo 



If the Percus-Yevick direct correlation function |33[ was used for cq one would end up with 
the Percus-Yevick result for the free energy in the special case of hard spheres (for which 
obviously a{oj) = l/(47r)). Instead, in accordance with Lee p4| , ^5| , we choose 



A <"> - S(T^ < 10 > 



which follows from the requirement that Eq. (|8|) yields the Carnahan-Starling expression |36 
for the free energy of hard spheres which is known to be more accurate at high densities. 
Note that fo{v) does not depend on the shape of the particles which enters Eq. (|^) only via 
the distance a(u 12 , uj, uj') of closest approach. 

In view of the molecular symmetry the orientational distribution a(u) is expected to be 
axially symmetric so that 

oo 

27ra(c<j) = a(cos 9) = aiPi(cos 9); (11) 

1=0 

Pi(x) are the Legendre polynomials. The excluded volume for fixed orientations u and uj' is 
given by 

v exd (u,uj') = w ex . c/ (cos7) = - / du 12 <J 3 (u 12 , uj,u') = 2_.viPi(cosj) (12) 

6 J 1=0 

where 7 denotes the angle between the directions uj and uj'. Insertion of the expansions 
Eqs. flllD and ([ID into Eq. © leads to 



^ = ^/.WE(^) v,al (13) 

As mentioned above this expression reduces to the Carnahan-Starling formula in the case 
of hard spheres. On the other hand in the limit 77 — ^ one recovers the first two terms of 
the virial expansion used by Onsager This limit is especially helpful for very elongated 
particles for which the isotropic-nematic transition occurs at very low packing fractions. 

B. Dipolar interaction 

The dipolar contribution to the free energy is treated in the so-called modified mean-field 
approximation f37H39f . 



Fdip = ~2/3 J d * rdudZr ' duJ ' a ( UJ ) a ( UJ ')®( r ^ ~ ^(^i2,^,^'))fdi P (r 12 ,uj,Lj'), (14) 

with the Mayer function fdi P = exp(—f3wdi P ) — 1 which is cut off at contact through the 
Heaviside function 0(r 12 — er(^ 12 , u, 00')). This expression follows from using the low- 



density approximation g rj exp(— f3w) for the pair distribution function. The Mayer function 
can be expanded as 



1 f -/3m 2 \ n 

fdip(ru,u,uj') = ( - a ) *'Vi2,w,^')- (15) 

n=l ' ^ 



ri2 3 



Due to the slow decay as function of r 12 the first term in this series requires particular 
attention. For this so-called /ong- ranged term 



F dt R) = ~P 2 / dudw'a{u)a{u') / rf 3 rdV6(r 12 - a{u 12 ,u,u'))w dip {r 12 ,u,u') (16) 



1 

2' 

the spatial integrations have to be analyzed carefully by first considering a fluid confined to 
a finite volume V and then performing the thermodynamic limit for a fixed shape of this 



volume. As has been shown in Ref. [20] the result does depend on the shape of the sample. 



For an ellipsoidal volume of aspect ratio k it was found that 
F (LR) 8vr 

— — = — p 2 m 2 a\(D{k) — 1/3) (for spherical particles) (17) 
V 9 



where D(k) is the demagnetization factor (see Eqs. (3.22) and (3.24) in Ref. 0]). For 
nonspherical particles the spatial integrations can be separated into contributions with < 
R c and r 12 > R c where R c is a fixed distance larger than the maximum of cr(uJi2, ui, to'). Since 
the result in Eq. QT7D does not depend on the particle size it can be adopted for the latter 
contribution (7*12 > R c )- The remaining integral (r±2 < R c ) can easily be evaluated since in 
this case the kernel is effectively short-ranged as function of r 12 yielding 

dip - p 2 m 2 al(D(k)- 1/3) + ^- / du;du/a(u;)a(u/)^(cos 7 ) (18) 



V 9 r iV v ' ' ' 2(3 
with 

g(^)(cos 7 ) =(3m 2 / duj 12 w(uj 12 ,uj,uj')\n Uc = V g { L ^P t (co8 7 ). (19) 

(As it should be q( LR ^ does not depend on the arbitrary parameter R c due to 
j duji2w(uji2-, uj, lo 1 ) = 0.) As mentioned above here we consider only the limit of an infinitely 
long and thin sample (k — ► 00) for which D(k) = 0. For other sample shapes the equilibrium 
configuration of a ferromagnetic fluid, i.e. with a% 7^ 0, exhibits an inhomogeneous structure 
with a spatially varying axis of preferential orientation, analogouus to the domain formation 
in solid ferromagnets. It was shown in Ref. |21[ that if this effect is taken into account the 
free energy is independent of the sample shape (see also Ref. PH ) and has the same value 
as for a spatially homogeneous liquid in the limit k — > 00 and thus the same phase diagram. 
The other terms in Eq. fllil) , i.e. the short-ranged contributions, can be written in the form 



rp(SR) 2 

1 dip _ P 

V ~ 2/3 

with 



dudu'a(u;)a(u')q {SR \cos-?) (20) 



? (Sfl) (cos 7 ) = / d^f; < 1)n ^, 3 „jf m2)W ,^ > 12 ,a;,a/). (21) 
J 3(n - l)n! a in ^{uj^.uj.uj') 



With the definition 



qi = / d(cos 7 )^(cos 7 ) (g (ifl) (cos 7 ) + <^(cos 7 )) - 5^[3m 2 (22) 



-1 



(D(k = oo) = 0) one finally has 



V 



1 2 
2 P ^U + l 



(23) 



C. Total free energy and phase coexistence 



Since Eq. (|23|) has the same form as Eq. (jl3|) we can write the total free energy as in 
Ref. @: 



£ = £| 



/l oo 
dxa(x)ln(2a(x)) + p 2 V 
1 Z=0 

with (now density- dependent) coefficients 

pui = 



(24) 



2TTl) + 



(25) 



If the summation over I is truncated at I — L (in practice it turned out that L = 4 is sufficient 
to yield reliable results) the minimization with respect to the orientational distribution leads 
to 

L 



MX) 



") ~ exp y-pP ^2 UiO>iPi(x) 

so that (see Eq. (5.6) in Ref. @) 

21 + 1 5-\ d> xP l( x ) ex P ^]f=i( 2i + l)uiOiPi{ 



(26) 
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(27) 



2 exp \ -p(3 Yli=i( 2i + l)uiOiPi\ 

The phase diagrams follow from requiring the equality of the chemical potentials and the 
pressures at the coexisting densities p\ and p 2 - 

dF 



dF 
dp 



dp 



(28) 



p 2 ,a( 2 ) (x) 



dF 



F(p u a^(x))- Pl — 



dF 



= F{p 2 ,a ( - 2 \x))-p 2 ^- 



p 2 ,a( 2 H x ) 



The functions a^'(x) denote the corresponding equilibrium orientational distributions ob- 
tained from Eqs. ( f26|) and (p7|) . As discussed above three kinds of phases are considered: 



S 



isotropic liquid (or gas) (a(x) = 1/2), nematic liquid (ot(x) = a(—x), i.e. a\ = for odd I), 
and ferromagnetic liquid (ai ^ for all I). (The latter phase could also be called ferromag- 
netic nematic, but we do not use this phrase in order to avoid confusion.) 

The determination of the phase boundaries for second order phase transitions is presented 
in the Appendix. Note that due to symmetry reasons there can be no truly second order 



isotropic- nematic transition |41| . 

We now discuss the calculation of the coefficients vi and q\ for the two types of particles 
we are interested in, i.e., hard spherocylinders and hard ellipsoids. A spherocylinder consists 
of a cylinder of length L and diameter D with two hemispherical caps of the same diameter 
D. The excluded volume is given by (see, e.g., Ref. [p2] ) 

All 

Uexrf(cos7) = 2L>L 2 sin7 + 2itD 2 L H D 3 (29) 

3 

which leads to 

t; = -DL 2 + 2ttD 2 L + — D 3 , v 2 = ——DL, V4 = -—DL 2 . (30) 

(Higher order terms are neglected.) The odd coefficients vanish due to the presence of the 
symmetry plane of the (nonpolar) particles. In order to determine the coefficients q\ one 
needs the function a{u\2, to'). It can be inferred from the observation that the surface of a 
spherocylinder is the set of all points with distance D/2 from the line segment connecting the 
centers of the caps. Thus when two spherocylinders are in contact these line segments always 
have the distance D from each other. It is difficult to give a closed formula for a(u)i2,u),u)') 
because several cases have to be distinguished. Nonetheless a numerical code for the cal- 
culation of this function can easily be implemented. An algorithm for the closely related 
problem of the distance between two spherocylinders with given positions and orientations 
has been derived explicitly by Allen et al. [[|. 

For hard ellipsoids with two equal axes of length a± and one axis of length an Berne and 



Pechukas [S3] introduced the often used approximation [31 



a(ui2,uj,uj') = a ± 



x- 



cos 2 9 + cos 2 6' — 2x cos 6 cos 6' cos 7 



1 — x 2 cos 2 7 



-1/2 



(31) 
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where 7, 8, and 9' are the angles between the directions u and lu', uj and CJ12, and Ct/ and 
U\2, respectively, and 



x = 4—4 = 4-T- ( 32 ) 



k = cr\\/a±. This approximation is obtained by considering the overlap of two ellipsoidal 
Gauss distributions. From Eq. (|3lD one finds 



Att 9 /T — cos^ 7^\ 

U«Brf(cOB7) = -g-0-||(7.L ( ) • ( 33 ) 



These formulas enable one to determine the coefficients vi analytically and qi numerically. 
We have truncated the sum in Eq. (|2l|) at n = 30 and found no significant changes upon 
including further terms for all considered values of the parameters /3m, a\\/a±, and D/L. 



III. DISCUSSION OF THE PHASE DIAGRAMS 

A. Non-polar ellipsoids 

Figure [I] displays the phase diagram for nonpolar hard ellipsoids with aspect ratio k = 
<j\\/a±. The solid lines denote the coexistence densities of the isotropic and the nematic 
fluid as determined from the theory presented in Sec. |T[ They are in fair agreement with 
the results of Monte Carlo simulations (squares) ||. The two-phase region is always very 
narrow, as in real nematic liquid crystals. The triangles indicate the liquid-solid transition 
found in the simulations, which cannot be described by the present theory. (The description 
of the solid phase requires a weighted-density-functional theory |j43| .) For aspect ratios near 
k — 1 the isotropic-nematic transition is preempted by freezing. Within the approximation 
in Eq. ( pi]) the physical properties are invariant under the transformation k — > 1/k. This 
behavior is satisfied very well by the simulation results, too. 
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B. Polar ellipsoids 



In this subsection we discuss the case that the hard ellipsoids are endowed with a point 
dipole of strength m oriented along the symmetry axis of length an . In the following we use 
the dimensionless reduced temperature T* = kBTa\/m 2 and the volume fraction 77 = pv^ 



transition is shifted to lower densities upon lowering the temperature. Below T* pa 0.2 the 
density gap increases dramatically and evolves into a broad coexistence region between an 
isotropic gas and a nematic liquid. Thus in this fluid gas-liquid coexistence is not terminated 
by a critical point. For rj > 0.265 lowering the temperature leads to two phase transitions 
(isotropic fluid— > nematic fluid and nematic fluid^gas-liquid) whereas for rj < 0.265 there is 
only the gas-liquid transition. As shown in Fig. |3] for k = 2 a ferromagnetic liquid appears 
in the medium temperature range. This phase turns into a purely nematic phase along the 
dotted lines of critical points. Whether the disappearence of the ferromagnetic order at low 
temperatures is an artefact of the approximations or not needs to be checked by alternative 
techniques. The behavior of the order parameters ol\ along different thermodynamic paths is 
displayed in Figs. |] and |5|. Figure |] illustrates their density dependence along two isotherms. 
For T* = 1.3 the ferromagnetic order parameter a,\ vanishes at the critical density with a 
square root singularity in accordance with the presently used mean-field theory while the 
nematic order parameter a<i exhibits a small break of slope not visible on the scale of the 
figure. At a lower density vanishes discontinuously at the first-order nematic-isotropic 
transition. For the lower temperature T* = 0.6 the ferromagnetic phase transforms directly 
into the isotropic phase. Figure [| shows the order parameters as function of T* at the fixed 
density rj = 0.75 demonstrating the loss of the ferromagnetic order at low temperatures. 
In the intermediate ferromagnetic region the orientational distribution a(cos6) typically 
exhibits two maxima at 6 = and 9 = n with the higher one determining the sign of 
the spontaneous magnetization. Their heights become equal at the nematic-ferromagnetic 




Fig. [I]) corresponds to the limit T* — > 00. For k = 3 (Fig. |2|) the first-order isotropic-nematic 
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transition. In contrast to Baus and Colot |24| , who use a different density-functional theory, 
we do not observe a phase transition between phases with one and two maxima in the 
orientational distribution. 

Upon further lowering the aspect ratio (Fig. ^ another qualitative change of the phase 
diagram occurs: for intermediate densities the isotropic-ferromagnetic transition becomes 
continuous. Two tricritical points arise where the character of this transition changes from 
first to second order. We remark that for high temperatures the ordered phases occur only 
at such high densities that they will certainly be preempted by a solid phase which is not 
captured by the present form of the density functional. 

Finally in Fig. ^ we present the phase diagram for dipolar hard spheres obtained from 
the present density-functional theory. There is no gas-liquid transition between isotropic 
fluids, but a coexistence of an isotropic gas with a ferromagnetic liquid which at a tricritical 
point (T£,r)t) = (0.600,0.198) changes into a continuous transition. As for the Stockmayer 
fluid |2(| the stability of the ferromagnetic phase is considerably overestimated as compared 
with simulations ]T^,T3] in which the isotropic-ferromagnetic transition has been detected 



at T* = 0.16 and rj « 0.4. 

For oblate particles (/t < 1) we obtain a similar series of phase diagrams (Figs. |^ and [|) 
as for k > 1 but without the loss of ferromagnetic order and the reentrance of the nematic 
phase at low temperatures. The comparison between Figs. || and ^ and between Figs. |9| 
and |] shows that the formation of the ferromagnetic phase is favored by oblate particles as 
compared with elongated ones. This is confirmed in Fig. [H^ which displays the phase diagram 
in the (/«, 77) -plane for a fixed value of the squared dipole moment per ksT normalized to 
the particle volume: [3m 2 /v^ = 6/n. The isotropic-ferromagnetic critical density increases 
with increasing elongations. An intuitive explanation for this observation is provided by 
an examination of the lowest energy configuration of two dipolar particles. The interaction 
energy at contact of a nose-to-tail arrangement is — 2m 2 /(K,a±) 3 while that of an antiparallel 
side- by-side arrangement is —m 2 ja\. Obviously the former configuration becomes more 
favorable as compared with the latter one if k is decreased implying a stronger tendency 
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for long-ranged ferromagnetic order. For large elongations the orientationally ordered phase 
becomes nematic while this does not happen for very oblate particles, for which instead 
a gas-ferromagnetic coexistence occurs at the temperature considered. Thus the k <-* 1/k 
symmetry of Fig. [TJ is lost in Fig. [TO as a result of the additional dipolar forces. 



C. Dipolar spherocylinders 

For dipolar spherocylinders one finds the same series of phase diagrams as function of 
the aspect ratio L/D as for elongated ellipsoids. An example is given in Fig. |ll| (compare 
Fig. As for the nonpolar ellipsoids the location of the isotropic-nematic transition for 
m = agrees very well with corresponding simulation results ||[8| which are available only 
for L/D = 5. 



D. Comparison with previous results 

In the following we compare our results with other theoretical investigations and numeri- 
cal simulations of dipolar hard core particles. Baus and Colot [|23J have studied this problem 
with a different density-functional ansatz which utilizes the analytically known correlation 
function of dipolar hard spheres in the mean-spherical approximation. At high temperatures 
they also find the sequence of isotropic, nematic, and ferromagnetic phases in accordance 
with our results, followed by a transition to a second ferromagnetic phase (see the end of the 
first paragraph in Subsec. [II FBI) . However, they claim that an increase of the aspect ratio 



shifts the ferromagnetic phase to lower densities and higher temperatures in contrast to our 
findings and to the qualitative argument presented at the end of Subsec. |III B| . 

The familiar second order virial expansion of Onsager ffl has been applied to dipolar 



spherocylinders by Vanakaras and Photinos p2|. It is well known that this approximation 



deteriorates for decreasing aspect ratios which is the reason that for L/D = 5 and 



m = these authors find the isotropic-nematic transition at densities above the density of 
closest packing. In addition they do not treat the long-ranged dipolar interactions correctly 
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which leads them to the wrong conclusion that a ferromagnetic phase cannot be stable for 
any dipole moment or aspect ratio. They agree with our results in that an increase of m or 
a decrease of T lowers the isotropic-nematic coexistence densities for point dipoles located 
in the center of the hard particles. However, as they point out, this trend is reversed if the 
point dipoles are located off center near the ends of the particles. 

The correlation functions obtained from the hypernetted-chain integral equation theory 
have been used by Perera and Patey [^3j for the investigation of dipolar hard ellipsoids. For 
k = 3 they find transition densities which are about 20% lower than the ones presented here 
and the same behavior of these densities upon increasing the dipolar strength. However, 
for k = 1/3 these authors were able to obtain a solution of the integral equation only for 
rather small dipole moments. They, too, did not take into account the correct long-ranged 
contribution to the free energy and did not detect a ferromagnetic phase for the aspect 
ratios under consideration. However, in a subsequent publication |18| they state that this 
conclusion remains unaltered if this fault is corrected. 



Vega and Lago have used a density-functional theory similar to the present one, 
but they incorporate a more sophisticated equation of state for the nonpolar isotropic fluid 
and treat the dipolar interactions within perturbation theory. The same remark as above 
concerning the treatment of the long-ranged interaction applies also to their work. Their 
results for the isotropic-nematic coexistence densities of dipolar hard sphero cylinders with 
aspect ratio L/D = 5 are in fair agreement with our findings for /3m 2 /D 3 between and 
4. However, they were unable to find a solution for the orientational distribution function 
for /3m 2 /D 3 = 6 which is close to the parameter range where within our approximation the 
isotropic-nematic transition broadens into a gas-nematic coexistence. 



Weis, Levesque, and Zarragoicoechea |27],[28|,|14]] have performed Monte Carlo simulations 



of dipolar spherocylinders with L/D = 5 for centered and off-centered dipole moments 
forming different angles with the cylinder axis. These authors were mainly interested in 
the structure of the smectic phase and did not determine the isotropic-nematic transition 
point. For {(3m 2 / D 3 ) 1 / 2 = 2.449 and a longitudinal centered dipole moment they report |27 
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the occurrence of an isotropic state at rj — 0.356 and a nematic state at r\ — 0.441 which 
is in accordance with the coexistence densities rji = 0.356 and i]^ = 0.375 obtained from 
the present theory. These authors have also studied ellipsoids with k = 3 [T2^] and found 
no isotropic-nematic transition for values of (/3m 2 /crj_) 1//2 in the range between and 3, in 
contrast to our results shown in Fig. |2|. For m = this finding is also at variance with Frenkel 
and Mulder's results ||, who did observe an isotropic-nematic transition with rjj = 0.507 
and i]n = 0.517. We conclude that further simulations of the isotropic-nematic or possibly 
the isotropic-ferromagnetic transition for different aspect ratios would certainly be helpful. 



Recently McGrother and Jackson [^4[ have published an extensive Monte Carlo study 
of the liquid-vapor coexistence in a dipolar hard sphero cylinder fluid. In agreement with 



previous work [45,17 they find that there is no liquid- vapor coexistence in the spherical 
limit L/D — > as well as for very elongated particles L/D 3> 1 due the formation of 
chains of nose-to-tail (small L/D) or side-by-side (large L/D) aggregated particles (see also 
Ref. |46|]). Only in an intermediate range around Lj D = 0.25 phase separation into isotropic 
gas and liquid does occur. This chain formation is not captured correctly by the present 
theory which predicts gas-ferromagnetic or gas-nematic coexistence for all values of L/D. 



IV. SUMMARY 

By applying density-functional theory for the description of orientational order in dipolar 
fluids consisting of hard nonspherical particles we have obtained the following main results: 

1. The location of the isotropic-nematic transition of nonpolar hard particles as a function 
of the aspect ratio has been obtained in good agreement with other theories and 
simulation results (Fig. [1]). 

2. The addition of a longitudinal point dipole at the centers of the particles induces a 
decrease of the coexisting densities of the isotropic-nematic transition and it leads to 
gas-nematic coexistence for large values of the dipole moment or at low temperatures 
(Fig. §. 
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3. There is also a ferromagnetic liquid phase provided the particles are not too elongated 
and the dipole moment is sufficiently strong. This phase is reached from the nematic 
or isotropic states by continuous or weakly first-order transitions, depending on the 
temperature and the particle aspect ratio. 

4. In accordance with the mean-field character of density-functional theory, at the con- 
tinuous phase transitions (Curie points) the magnetization vanishes according to a 
square root power law and the nematic order parameter exhibits a small break of 
slope (Figs. £| and §). 

5. Within the present theory the phase diagram of dipolar hard spheres (Fig. [7]) comprises 
an isotropic and a ferromagnetic fluid with first (second) order phase transitions at 
temperatures below (above) a tricritical point. 

6. Oblate particles exhibit a similar phase behavior as elongated ones, but a stronger 
tendency for the formation of the ferromagnetic phase. In contrast to the nonpolar 
case the phase diagrams are not approximately symmetric with respect to the trans- 
formation k «-> 1/ac, where k is the aspect ratio of the uniaxial ellipsoidal particles 
(Fig. 0). 

7. For dipolar hard spherocylinders we find an analogous series of phase diagrams as for 
elongated ellipsoids. 



APPENDIX: CRITICAL DENSITIES 

We rewrite the contributions to those free energy in Eq. fl24|) which depend on the 
orientational distribution as (here and in the following all integrals over x and x' are to be 
taken over the interval [—1,1]) 

— I dxa(x) ln(2a(x)) H — p I dxdx' a(x)d(x > )K(x, x') (Al) 



V (5 J w v v » 2 ' 
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with 



K(x,x') = ^t^ u^P^x'). (A2) 



1=1 



Following van Roij et al. |47j we use a kind of bifurcation analysis in order to determine the 
critical density of the nematic-ferromagnetic or isotropic-ferromagnetic phase transition for 
a given temperature. The minimization of Eq. ( |A1| ) yields 



ln(2a(x)) +0p dx'a(x')K(x,x') = v (A3) 



where the Lagrange multiplier v is determined by the normalization j dxd(x) = 1. One con- 
siders a small ferromagnetic perturbation a>i(x) of a solution a>o(x) with nematic symmetry, 
i.e. a (x) = a (— x). If the expansions a(x) = a (x) + edi(x) + . . . and v — i/ + ev\ + . . . 
are inserted into Eq. (|A3|) the term linear in e gives 

^rl + [3p [ dx' aAx'^Kixix') = v x . (A4) 
a (x) J 

Integrating Eq. ( |A4|) over x and using the relation J dxK(x,x') = yields 
a x (x) 1 f ,di(x') 



a Q (x 



- f dx'^4 + (3p [ dx 1 aAx^Kix.x') = 0. (A5) 
2 J a (x') J 



Since both a (x) and de(x) must be of the form exp(^2f =0 "fiPi(x)) (see Eq. (EH)), with the 
nematic solution do(x) containing only even indices I, we make the following ansatz for the 
small perturbation: 

L/2 

6li(x) = d (x) ^ 72n-l-P2n-l(^) (A6) 
n=l 

where terms of the order 72n-i72m-i and higher due to the exponential form (see above) 
have been neglected. Due to Eq. (|A5|) the coefficients 7; with odd I satisfy the equation 

L/2 
l'=l 

with 

(21 + l) 2 r 

A w = K 2 J Ppui / dxd^P^Pvix). (A8) 
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Both the nematic solution cto(x) and the L/2 x L/2 matrix A depend on the density. The 
critical density is reached if one of the eigenvalues of A equals —1 giving rise to a nontrivial 
solution of Eq. (|A7|). For L = 2 this condition reduces to 

9 f 

-(3pui / dxao(x)x 2 = —1. (A9) 



2 

For L = 4 one obtains after some algebra 



1 + A u + A 33 + A U A 33 - A 13 A 31 = 0. (A10) 



From these equations together with the numerical solution of Eq. (|27|) the nematic- 
ferromagnetic critical density can be determined. 

For the isotropic-ferromagnetic transition the unperturbed solution is ao(x) = 1/2 so 



that Eq. (|A9|) reduces to the known result (see Eq. (7.10) in Ref. [20]) 



Ippui = -1. (All) 

In this case the matrix A is diagonal for general L: Au> = (21 + l)/2(3pui5^i>. Thus the above 
bifurcation condition yields 

^Y^Pfmi = -I- (A12) 

The actual values of the coefficients ui are such that the lowest density for which Eq. (|A12|) 
is satisfied always corresponds to / = 1, so that Eq. ( |A11| ) is valid for all L. 
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FIGURES 

FIG. 1. Phase diagram of nonpolar hard ellipsoids with two equal axes of length a± and one axis 
of length ail and an aspect ratio k = cr^/aj_. r\ = pv^ ' is the volume fraction where p is the number 
density of the particles and v<® is their individual volume. The lines denote the coexisting volume 
fractions at the first-order isotropic-nematic transition as obtained from density-functional theory. 
The squares and triangels are simulation results for isotropic-nematic and liquid-solid coexistence, 
respectively. For the former there is good agreement with density-functional theory. The dotted 
lines are guides to the eye. The isotropic-nematic transition is accompanied by only a small density 
discontinuity which decreases for k — > 1. The present density- functional theory is not suited to 
describe the freezing transition. The simulation data support the symmetry k «-> 1/k discussed in 
the text. 

FIG. 2. Phase diagram of dipolar hard ellipsoids with aspect ratio k = 3 in the tempera- 
ture-density plane. The weakly first-order isotropic-nematic transition at high temperatures broad- 
ens into gas-nematic coexistence at low temperatures. Here and in the following figures two-phase 
coexistence regions are shaded. 

FIG. 3. The same as in Fig. ^ but for k = 2. Besides the isotropic (I) and nematic (N) phases 
a ferromagnetically ordered liquid (F) occurs in the intermediate temperature range. The dotted 
and solid lines denote second- and first-order transitions, respectively. The density gap of the I-N 
transition at rj ~ 0.68 cannot be resolved on the present scale. The continuous F-N transition at 
low temperatures intersects the first-order transitions at a critical end point. The inset shows that 
in a narrow temperature range the high temperature continuous N-F transition is turned into a 
weakly first-order transition generating a tricritical point and an I- N-F triple point. The two-phase 
coexistence regions are shaded. For reasons of clarity this shading is omitted for I-F coexistence. 



FIG. 4. The first two orientational order parameters a\ (see Eq. (11)) along two isotherms in 
Fig. |3] for k = 2. The gaps between the black dots indicate two-phase regions. The dotted vertical 
lines indicate discontinuities. 
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FIG. 5. Temperature dependence of the first four orientational order parameters ai along a 
thermodynamic path of fixed density in Fig. || for k = 2. The order parameters with odd indices I 
vanish in the low- and high-temperature nematic phases but are nonzero between the lower critical 
point T*i = 0.226 and the upper critical point T* 2 = 1.34 corresponding to ferromagnetic order, ai 
and ct!4 are nonzero for all temperatures and decrease for increasing temperature. At T*y and T* 2 
they exhibit a break of slope. 

FIG. 6. Phase diagram of dipolar hard ellipsoids for k = 1.5. The line styles have the same 
meaning as in Fig. ||. In the high temperature region the gap between the two coexisting densities 
cannot be resolved on the scale of the figure so that only a single solid line is visible. Between two 
tricritical points the isotropic-ferromagnetic transition is continuous. The two-phase coexistence 
region at low temperatures is shaded. 

FIG. 7. Phase diagram of dipolar hard spheres. Within the present approximation only one 
isotropic fluid and a ferromagnetic liquid are stable at any temperature. Below (above) the tricritical 
point the phase transition is discontinuous (continuous). 

FIG. 8. Phase diagram of oblate dipolar ellipsoids with k = 2/3. The meaning of the solid and 
dotted lines is the same as in Fig. |||. In contrast to Fig. [6| there is no reentrant nematic phase at 
low temperatures. 

FIG. 9. Phase diagram for oblate dipolar ellipsoids with k = 1/3. See also the caption to Fig. |3[ 

FIG. 10. Phase diagram of hard ellipsoids for a fixed reduced dipole moment in units of ksT 
and the particle volume v^°\ The stability of the ferromagnetic phase decreases with increasing 
aspect ratio k until it finally turns into a nematic phase at high densities and into an isotropic 
fluid at low densities. Coexistence of an isotropic gas and the ferromagnetic liquid occurs for oblate 
particles. There is no k <-> 1/k symmetry as in Fig. |l| due to the presence of the dipolar interactions. 
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FIG. 11. The phase diagram for polar spherocylinders with an aspect ratio D/L = 1 exhibits 
a similar behavior as that for polar ellipsoids with k = 2 (see Fig. ||). 
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Fig. 1 
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Fig. 3 
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Fig. 4 
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Fig. 5 
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Fig. 6 
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